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I. 



INTRODUCTION 



The history of digital computing is characterized by 
the near exponential progress in hardware, operating sys- 
tems, programming languages, and numerical computing algor- 
ithms. In contrast, the science of pseudo-random number 
generation has virtually stood still for the last twenty 
five years in spite of literally hundreds of publications 
on the subject. 

D. H. Lehmer proposed the linear congruential scheme 
for the generation of pseudo-random numbers on a digital 
computer in 1951. This scheme has been the singularly most 
popular method ever since. Recently, feedback shift-register 
methods have been proposed for pseudo-random number generation 
but they have not achieved as widespread use. 

It is essential to recognize the fact that neither method 
produces truly "random" numbers by any definition of random- 
ness. Consequently, a great deal of effort has been invested 
in examining the sequences produced by random number genera- 
tors to establish whether they meet statistical tests for 
randomness for relatively short samples. Given that a 
generator has "passed" an acceptable battery of statistical 
tests, its sequence is considered adequate for use in simulation 
experiments . 

An elaborate collection of statistical tests has been 
developed over the years and Knuth [Ref. 8] provides an 
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excellent discussion of these tests. They are, however, not 
sufficient for certifying the quality of a generator. The 
report by Learmonth and Lewis [Ref. 9] attests to this fact 
in that for certain generators known through use to be poor, 
the statistical tests were not conclusive in identifying 
their weaknesses. 

The recent work on testing Lehmer congruential random 
number generators has concentrated on characterizing the 
entire periodic sequence. The two tests proposed, which are 
essentially similar, are the lattice test and the spectral 
test. These tests provide deterministic rather than statis- 
tical criteria for rating the quality of the sequences. 

The thesis to be examined here is that although these 
two tests are a significant advance in the testing of Lehmer 
congruential random number generators, they are formulated 
incorrectly and the computational algorithms to implement 
the tests are unnecessarily complicated. 

It will be shown that the proper space in which to 
characterize finite periodic sequences is an algebraic group 
rather than the space of the natural numbers or the real 
line as has been previously assumed. The development to 
follow will concentrate on the class of primitive root/prime 
modulus random number generators which form a more general 
algebraic system, a finite field. 

The special properties of finite field arithmetic will 
be developed and several new theorems will be presented 
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regarding primitive roots of prime fields. Using these 
properties of finite fields, the structure of linear con- 
gruential sequences on finite fields will be examined. 

The conventional lattice and spectral tests will be 
defined and their algorithms outlined. These two tests 
will then be recast for finite field assumptions. New 
algorithms will be sketched as potential replacements for 
the lattice and spectral tests. 

While the surface will only be scratched, the underlying 
theory presented will form the basis for future research 
in random number generator testing. 
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II. FINITE FIELDS 



In this section, it is intended to make a critical 
reevaluation of the assumptions underlying the mathematical 
characterizations of linear congruential sequences. Both 
the spectral test of Coveyou and Macpherson and Knuth and 
the lattice test of Marsaglia and Ahrens and Dieter makes 
the very important assumption of an infinite sequence of 
integers. As will be shown, for the spectral test this 
assumption is critical for the development of the equidis- 
tribution properties of linear congruential sequences. 
Although the lattice test does not require this assumption, 
it is inherent in the development of the computational 
algorithm. 

Abandoning this assumption and using the more appro- 
priate assumption of a finite sequence, it will be seen that 
both the theoretical as well as computational development 
will be considerably more clear and concise. 

A. FINITE FIELDS 

The integers resulting from the linear congruential 
equation defining a random number generator from an alge- 
braic group. A group is defined as a set of elements upon 
which is defined a binary operation called multiplication 
with the following properties: 

1. the operation is closed, that is if a and b are 
elements of the group then c = ab is also in the group; 
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2. the operation is associative, that is 



a (be) = (ab) c ; 

3. there is an identity element in the group, usually 
denoted by the integer 1 such that la = a holds for every 
element a of the group; 

4. for each element a of the group there is an inverse 
denoted a ^ such that a ^ a = 1 . 

The residues which are relatively prime to a general 
modulus m under the operation of multiplication modulo m 
form a group as defined above. When the group operation is 
multiplication, the group is called a multiplicative group. 
Similarly, an additive group may be defined with the group 
operation being addition. The inverses in the multiplicative 
group fill the role of division. For additive groups, the 
inverse is subtraction and the identity element is 0, i.e., 
a + 0 = a. In either case, if the group operation is also 
commutative, ab = ba or a + b = b + a, the group is given 
another adjective. Abelian. 

To recapitulate, for a random number generator 
X^ + ^ E a X^ (mod m) where m is composite, the integer residues 
form a finite Abelian multiplicative group under the opera- 
tion of modulo multiplication. For the present case the 
residues of the random number generator X^ + ^ E a X^ (mod p) , 
where p is prime, form a more general algebraic system, 
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specifically they form a finite field. A field is a gener- 
alization of a group wherein the group operations of addition, 
subtraction, multiplication, and division are all defined 
either as group operations or their respective inverses. 

It is the arithmetic properties of finite fields that 
will reveal the structure of the linear congruential sequences. 

B. FINITE FIELD ARITHMETIC 

Given a finite field as the algebraic structure on which 
the residues from the random number generator are defined, 
it remains to examine how the ordinary arithmetic operations 
work in this system. 

Modulo addition works as would be expected. In ordinary 
addition over the integers, the only solution to the equation 
a + x = 0 is x = -a. In finite field arithmetic over the 
positive integers the solution is given by another positive 
integer called the additive inverse. For example, in the 
finite field formed by the prime 17, the additive inverse 
of the integer 3 is not -3 but 14 

3 + 14 = 17 5 0 (mod 17) 

Therefore x = 14 satisfies a + x = 0 in this finite field. 

Modulo multiplication also works as in conventional 
arithmetic. The multiplicative inverse, which plays the role 
of a divisor is somewhat differently defined. In ordinary 
arithmetic over the integers the solution to ax = 1 is 
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uniquely x = 1/a = a In finite field arithmetic, the 

role of a ^ is taken by another positive integer in the 
field. Again referring to the field formed by the prime 
17, the solution to 3 x = 1 is x = 6, i.e., 3x6= 18 = 1 
(mod 17) . Therefore 3 ^ = 6 in this field. 

Although subtraction and division are defined as the 
inverses of the field operations of addition and multipli- 
cation, the algebra of equations over the finite field 
becomes more difficult. The absence of negatives makes 
solutions to simple equations such as 3 - x = 7 impossible. 

Since the labelling of the field elements as positive 
integers was arbitrary, it is possible to relabel the ele- 
ments more conveniently. One such relabelling, is as 
follows , 

{- |(p-l), - |(p-3) , ..., -1, 0, 1, ..., |(p-3) , y(p-l)} 

Here the elements previously labelled with integers greater 
than ~(p-l) are mapped into their additive inverses which 
are less than or equal to j(p-l) and of opposite algebraic 
sign. For instance, in the field formed by the prime 17, 
the positive integer 14 is mapped into the element -3. This 
mapping preserves the additive inverse of 3 since 
3 + (-3) = 0 as required. This new notation will be called 
the primitive mark notation. 
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C. PRIMITIVE ROOTS 



Before examining the structure of linear congruential 
sequences, it will be instructive to examine primitive roots 
of primes in the context of finite fields. As previously 
indicated, for the type of random number generators to be 
examined here, namely the primitive root/prime modulus type, 
the, multiplier is required to be a primitive root to guar- 
antee a full period. In finite field terminology primitive 
roots are termed generators of the field. 

The term generator derives from the finite analogy of a 
generating function. Any primitive root of the prime finite 
field may be used in defining the generating function on 
the finite field. The term primitive root also has the very 
special connotation of primitive root of unity. A primitive 

i 

root serves for a finite field the same special properties 
that the primitive root of unity exp {2irik} serves for the 
complex field. In fact, the (infinite) field of complex 
numbers is the only other field which possesses a primitive 
root of unity. 

While it is known that every prime possesses at least 
one primitive root, the only other significant fact that 
number theory offers is that each prime p contains precisely 
<j>(p-l) primitive roots. Here <j)(p-l) is Euler's totient 
function which is defined to be the number of integers less 
than and relatively prime to (p-1) . When constructing a 
primitive root/prime modulus random number generator the 
greatest problem lies in discovering a primitive root of the 
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given prime. The procedure usually followed is to use a 
trial and error technique on small integers and then to 
apply the definition of a primitive root to verify if the 
number is in fact primitive. This actually is not an easy 
task since showing that p-1 is the least positive exponent 
such that a^ ^ = 1 (mod p) requires considerable computation. 
Tables do exist listing several primitive roots for some 
primes and this is a great help. Knowledge of one primitive 
root provides a rather easy method for finding others. A 
simple, but potentially very lengthy algorithm for finding 
all the primitive roots of a given prime is as follows. 

Algorithm A : 

Al. Factor p-1 into its prime factors. (A sieve 
method may be used.) 

A2 . Select a known primitive root, say a? set k = 1. 

A3. Divide k by each of the prime divisors of p-1. 

If any of the prime divisors divides k evenly, 
i.e., no remainder, then a is not a primitive 
root, go to A4 . Otherwise print a as a 
primitive root. 

A4. If k is greater than or equal to p-2, stop, all 
primitive roots have been found and printed; 
otherwise k = k+1, go to A3. 

Viewing primitive roots as generators of the finite field 
of characteristic p, it is possible to use the properties of 
finite fields to discover further knowledge about other 
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primitive roots. All finite field elements may be represen- 
ted as power residues of a primitive root of the field, that 
is, every element may be represented as a (mod p) for a 
unique k. The multiplicative inverse of any finite field 
element is found easily using the following lemma. 

Lemma 1 : The multiplicative inverse of a finite field 

element (not necessarily a primitive root) , a (mod p) , is 
simply a^ ^ ^ (mod p) . 

„ - k (p-l)-k _ k (p-1) -k , , . 

Proof : a a ^ = a a ^ a (mod p) 

= 1 (mod p) 

since a^ ^ =1 (mod p) by definition of a primitive root. 
Q.E.D. 



This lemma leads to the following theorem concerning the 
multiplicative inverses of primitive roots. 

Theorem 1 : If a is a primitive root of the prime p, then 

its multiplicative inverse a ^ = a^ ^ (mod p) is also a 
primitive root. 

Proof: The multiplicative inverse of a is, by Lemma 1, 



a*l 5 atP - 1 ’- 1 



a <P-2> 



(mod p) . 



—*1 ( D 2 

To show that a E a vp ) (mod p) is a primitive root of p, 
it must be shown that (p-1) is the least positive exponent 
such that 
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(p~2) , (p-1) = 



1 (mod p) . 



(a' 1 ) ‘P- 1 ’ i (a 



Hence , 



(a -l, (p-l) 5 (a (p-2) , (p-1) E a (p-2)p a -(p-2) 

5 a (P'D (p-1) -1 a “ (p-2) 

= 1 a’ 1 a'<P- 2 > 

= a ( P— 2 ) a -(p-2) 

= 1 (mod p) . 

Assume (a^ = 1 (mod p) with 0 < q < (p-1) , then 



(a 



(p-2) 



-q = a (p-Dq a -q = a -q 



(mod p) 



Applying Lemma 1 again 

a ^ = a^ ^ ^ (mod p) . 

If a ^ is congruent to 1 modulo p this would imply a^ ^ ^ 
is also congruent to 1 and this contradicts the assumption 
that a is itself a primitive root, i.e., (p-1) is not the 

least positive exponent such that a^ ^ = 1 (mod p) . Q.E.D. 
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If the multiplicative inverse of a primitive root is also 
a primitive root then it would be logical to examine additive 
inverses. Unfortunately the case for additive inverses is 
somewhat more restrictive. To examine the additive inverses 
of primitive roots it is necessary to separate the primes 
into two classes, specifically the classes p = 1 (mod 4) 
and p = 3 (mod 4) for primes p > 2. It is trivial to see 
that these are the only two classes into which the primes 
fall modulo 4 . 

Theorem 2 : If a is a primitive root of the prime p > 2, 

then the additive inverse of a is also a primitive root of 
p if and only if p = 1 (mod 4) . 

Proof : In the unambiguous notation of primitive mark repre- 

sentation, the additive inverse of a primitive root a is 
-a. To show that -a is a primitive root of p when p = 1 (mod 4) , 



(p— 1) (mod 4) 



(-a) 



1 (mod p) . 



To show that -a is not a primitive root of p when p e 3 (mod 4) , 



, , (p-1) (mod 4) , ,2 2 . , , . 

(-a) r e (-a) E a El (mod p) . 



Since it was assumed that a was a primitive root of p, the 
above contradicts the primitivity of a, i.e., p-1 is not the 
least positive integer such that a^ ^ El (mod p) for 
primes p > 2. Q.E.D. 
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The two theorems just presented offer some insight into 
the distribution of primitive roots over a prime finite field. 
Theorem 2 is particularly instructive when considering the 
field integers in their primitive mark representation. Know- 
ing that the additive inverses are also primitive roots for 
certain primes, the <J>(p-l) primitive roots are distributed 
symmetrically about the origin at 0. Considering again the 
finite field formed by the prime (17 = 1 (mod 4)), there 
are 0(16) = 8 primitive roots distributed as follows: 

-8 -7 -6 -5 -4 -3 -2 -1 0 1 2 3 4 5 6 7 8 

The respective multiplicative inverses are indicated. 

For the case of primes p 5 3 (mod 4) , the situation is 
less appealing. Since the additive inverses are not primitive 
roots, the symmetry is lost. There is no particular pattern 
to the distribution of the primitive roots or their 
multiplicative inverses. 
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III. STRUCTURE OF LINEAR CONGRUENTIAL SEQUENCES 



The attempts to characterize the sequences emanating from 
linear congruential random number generators recently have 
focused on structural representations. The lattice test im- 
poses a grid on the generated n-tuples of points and charac- 
terizes the sequences in terms of the lengths of the sides of 
the smallest n-dimensional hypercube which can be constructed. 

The spectral test, which can be shown to be a variant of the 

\ 

lattice test, views the n-tuples as forming waves of hyper- 
planes and characterizes the sequence by the "frequency" of 
these waves . 

Marsaglia [Ref. 9] attempted to expose more of the repeti- 
tive structure of the generated sequences. Three new theorems 
will be presented here to further characterize the structure 
of these sequences. By exploiting the primitive mark notation 
for finite fields more insight will be gained. These new 
theorems will then be contrasted with Marsaglia' s work, lead- 
ing to a redevelopment of the lattice test and spectral 
test. 

A. THE FUNDAMENTAL SUBSEQUENCE 

It is known that the period of a primitive root/prime 
modulus random number benerator is (p-1) . Using the positive 
integer representation for the elements of the finite field, 
the sequence appears to have full period. The following 
theorem will establish the existence of a more fundamental 
half sequence. 
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Theorem 3 : The sequence of (p-1) integers generated by a 



primitive root/prime modulus random number generator, when 
expressed as primitive marks, consists of a fundamental sub- 
sequence of length -j(p-l) . This fundamental subsequence is 
followed by a replicate subsequence identical to the first 
except for opposite algebraic sign. These two subsequences, 
each of length i(p-l) constitute the effective period of 
length (p-1) . 

Proof : Without loss of generality, let the sequence begin 

with Xq = 1, then x^ = a, the primitive root multiplier. In 
th k 

general, the k element is x^ E a (mod p) . After generating 
■j(p-l) elements of this sequence, the next integer element 
is 

|(P-D 

x. =2 (mod p) . 

-j(p-D 

Now, by definition of a primitive root, a^ p =1 (mod p) , 
therefore 

x? = (a ?(P 1} ) 2 E a (p-1) = i (mod p) . 

i(p-D 

This implies that x. E ± 1 (mod p) . In fact, 

jtp-D 

x. = -1 (mod pf (the primitive mark representation of 

?<P-D 

(p-1) ) . Since it is known that each integer appears once and 

only once in the sequence and +1 appeared as Xq, therefore, 

beginning with element x-, , the same fundamental subsequence 

j(P-D 

of primitive marks repeats with opposite algebraic sign. Q.E.D. 
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For illustration, consider the following example using 



the prime modulus 17 and primitive root 3. 



Deviate ; 


X 0 


h X 2 


X 3 


X 4 


X 5 


X 6 


*7 


Positive 

integer 


1 


3 9 


10 


13 


5 


15 


11 


Primitive 

mark 


1 


3 -8 


-7. 


-4 


5 


-2 


-6 



X 8 *9 X 10 *11 *12 *13 ^4 ^5 ^6 

16 14 8 7 412 2 61 

-1 -3 8 7 4 -5 2 6 1 



FIGURE 3.1 X i+1 = a X ± (mod p) ; a = 3, p = 17 

For the lattice test, the consequences of this theorem 
are not particularly important. The grid formed by the lattice 
of n-tuples is very easily constructed using very few n-tuples 
and its regularity is apparent without having to generate the 
entire sequence. For the spectral test, however, the ffect 
would be a folding of the sequence onto itself. 

The next theorem relates the properties of multiplicative 
inverses to the resulting sequence generated by primitive 
root/prime modulus random number generators . 

Theorem 4 ; In the primitive root/prime modulus random number 
generator = a (mod p) , replacing the multiplier a by 

its multiplicative inverse, a \ results in a generated sequence 
which is precisely the reverse of the sequence generated by 

a . 
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Proof : The sequence produced by X^ + ^ E a X^ 

with X Q = 1 is X: = {1, a, a 2 , a (p_2) . 

The sequence produced using the inverse a 1 
X* = 1 is X*: = {1, (a -1 ) , (a -1 ) 2 , . .., (a -1 

(mod p) . By definition of a primitive root, 
therefore. 



(mod p) starting 
a ^ P D } (mod p) . 
starting with 

) (p- 2 ) / (a -vp-i) } 

a ( p -l) E 1 (mod p) 



1 = X* = X _ 1} E a 155 " 15 E 1 (mod p) . 



Similarly, 

a ^ = X£ = = a^ p 2 ^ E a^ p ^a ^ E a ^ (mod p) . 

By induction, for the general element X*, the corresponding 

element of X is X, , ,, , T 

(p-k-1) , (see Lemma 1) , 

-k v „ (p-k-1) _ (p-1) -k _ -k , , . 

a = X* = X^_ k _ 1 ^ = a ^ = a ^ a = a (mod p) . 

Q.E.D. 



Clearly, the sequence generated when the multiplicative in- 
verse is substituted, is precisely the reverse of the one 
generated using the original primitive root. Application of 
Theorem 3 will reveal that the fundamental subsequences are 
reversed also. 

Essentially, substituting the multiplicative inverse for 
the primitive root multiplier causes the random number generator 
to run "backwards." 
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For the special case of prime moduli such that 
p = 1 (mod 4) , substituting the additive inverse for the 
primitive r-ot multiplier results in further interesting 
properties of the sequence as the following theorem 
demonstrates . 

Theorem 5 : In the primitive root/prime modulus generator 

X^ + ^ = a (mod p) , where p = 1 (mod 4) , replacing a by 
its additive inverse, -a, results in a generated sequence 
with every odd numbered element of the sequence replaced by 
its additive inverse. Expressed in primitive mark notation, 
the fundamental subsequences are the same as for the multi- 
plier a except that the odd numbered elements have opposite 
algebraic sign. 

Proof : Without loss of generality, let the sequence 

X^ + ^ = a X^ (mod p) begin with Xq = 1, then the generated 
sequence is 



X: = {1, a, a 2 , ..., a^ p 2 ^ , a^ p ^ } (mod p) 



f(p-l) 

It was shown in the proof of Theorem 3 that a = -1 (mod p) , 

1 1 

hence a 2 ^ p ^ = a^ P ^ a = -1 a = -2 (mod p) . That is, 
a Z(P 1) +1 t j ie additive inverse of a. The sequence 
generated by X£ + ^ = (-a) X* (mod p) beginning with X* = 1 
is 



X*: = (1, (-a), (-a) 2 ,..., (-a)^ p 2 \ (-a) ^ P ~ 1 ^ } (mod p) . 
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The elements (-a) for k even is X* are obviously equal to 
k v 

the elements a in X. The elements (-a) for k odd, that is, 
the odd numbered elements of X* are the additive inverses 
of the elements a in X in primitive mark notation, 

+ (-a)^ = a^ - a^ = 0 (mod p) . 

The second part of the theorem follows directly from Theorem 
3. Q.E.D. 



Where Theorem 4 showed how to reverse the sequence using the 
multiplicative inverse, using the additive inverse effects 
a "half shuffle" of the sequence. 

B. MARSAGLIA 'S THEOREMS 

Marsaglia (Ref. 9] disucssed the idea of a "fundamental 
sequence" for the third generator X^ + ^ = a X^ + b (mod m) 
where m is composite. This theorem will be presented here 
and then contrasted with Theorem 3. 

Theorem 6 (Marsaglia [Ref. 9]) : The choice of b and the 

starting value Xg are of no consequence for the generator 

X^ + ^ = a X^ + b (mod m) , in the sense that if Xg , X^, X^ r 

... is any sequence generated by X^ + ^ = a X^ + b (mod m) and 

if Yq, Y^, Y^t ... is the fundamental sequence 0, 1, 1+a, 

2 

1+a+a , . . . , generated by Y^ + ^ = a Y^ + 1 (mod m) , then 

there are constants V and W such that X^ = V Y^ + W. In 
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Thus the sequence of the 



fact, V = Xg(a-l) + b and W = Xq . 
form X ±+1 i a X i + b (mod m) with X Q and b arbitrary (includ- 
ing b = 0) may be obtained from an affine transformation of 
the fundamental sequence. 

In addition, if the multiplier a is relatively prime to 
m, then the period of the sequence Xq, X^, X 2 , ... with 
= a X^ + b (mod m) is the period of the fundamental 
sequence for modulus m/d where d = gcd {m, X^ia-l) + b}. 

(End of theorem) . 

For the case of the primitive root/prime modulus genera- 
tor, m = p, a prime, a is a primitive root of p, and b = 0. 
With no loss of generality, for Xq = 1 

d = gcd {p, (a-1) } = 1 

since (a-1) and all other integers less than or equal to 
(p-1) are relatively prime to p. Therefore the period of 
this fundamental sequence is theperiod of the sequence for 
p/d = 1 or (p-1) as is already known. The sequence generated 
by X^ + ^ = a X^ (mod p) can be put simply into one-to-one 
correspondence with the sequence {4} in the theorem. In the 
case of the primitive root/prime modulus generator, Marsag- 
lia's theorem is completely unins true t ive , and by using 
positive integer notation, it actually misses the existence 
of the fundamental subsequences of Theorem 3 which uses 
primitive mark notation. 
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Marsaglia provides another theorem which attempts to 
uncover more detail about his fundamental subsequences. It 
is presented here. 

Theorem 7 (Marsaglia [Ref. 9]) : Let a be relatively prime 

to the modulus m and have multiplicative order 1 t. Then 
the fundamental sequence: 

Y g = 0, Y 1 = 1, Y 2 = 1+a, Y 3 = 1+a+a 2 , ... 

Y i+1 = 1 + a Y j_ (mod m) (1) 

is made up of a block {B} = {Yg, Y^, ..., °f t distinct 

residues of m, followed by translates of that block {B}, 
{B+c}, {B+2c}, ... where c = l+a+a 2 + ... +a^ fc 1 ^ (mod m) . 

The period of the sequence (1) is tr, where r is the additive 
order of c: r = m/gcd {c,m}. 

Applying this theorem to the primitive root/prime modu- 
lus case where m = p, a prime, b = 0, and a is a primitive 
root of p, it is clear that t = (p-1) by definition of a 
primitive root. Therefore, Marsaglia' s fundamental sequence 
is made up of a block {B} = {Yg, Y^, ..., 2)} °f the 

(p-1) distinct residues of p,i.e., the sequence constitutes 
one block of length (p-1) . For the additive order of c, 



integer 



Multiplicative order is defined as the least positive 
ger t such that a r = 1 (mod m) when gcd {m,a} = 1. 
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p- 2 , 

Z a J (mod p) 



j = 0 



1 - a'?' 1 ' 



(mod p) 



0 ; 



r = p/gcd {0,p} = p/p = 1 and tr = (p-1) . 

As in the case of Theorem 6, this theorem provides no 
information. Theorem 3 offers much more insight into the 
structure of the linear congruential sequence. The recog- 
nition of the finite field and using primitive mark notation 
provide the basis for the redevelopment of the lattice test 
and spectral test. 
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IV. THE LATTICE TEST 



The use of the lattice test is a fairly recent develop- 
ment even though its proposal as a means of characterizing 
linear congruential sequences goes back at least twelve 
years. Franklin [Ref. 6] first observed the hyperplane 
structure of n-tuples from a linear deterministic sequence. 

He was concerned, however, with infinite sequences. Jans- 
son's book [Ref. 7] reiterated the usefulness of the lattice 
and alluded to its possible use with finite, periodic se- 
quences . He did not provide an algorithm or any computational 
results . 

Marsaglia's famous paper [Ref. 10] focused widespread 
attention on the hyperplane structure of the sequences pro- 
duced by congruential generators. Papers by Beyer, Roof, 
and Williamson [Ref. 4] and Smith [Ref. 13] published the 
first algorithms and computational results. 

Recent results on the lattice structure of popular gen- 
erators come chiefly from the later paper of Marsaglia 
[Ref. 9] and the as yet unpublished book by Ahrens and 
Dieter [Ref. 2]. Both of these references contain compu- 
tational algorithms for performing the test. It is the 
difficulty of implementing these algorithms on a digital 
computer which inspired the work to follow. 
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A. THE ALGORITHMS 



The current algorithmic implementations of the lattice 
test will now be described with emphasis on the assumptions 
made. Marsaglia's algorithm, BEST2 , will be discussed first. 

In n-dimensional space, it is assumed that a set of 
points b^, , ..., b fi exist such that all points in the 

space are integral multiples of them, that is, 

{r^ + r 2 b 2 + ... + r n b n ; r^r^ ..., r n = 0,±1,±2, ...} 

constitutes the set of all points spanned by b^, b 2 , ..., b . 
These points then form a lattice in n-space. 

The points b^, b 2 , . . . , b n are actually n-vectors so that 
when viewed as rows of a matrix they form a basis for the 
points in n-space. Any set of n linearly independent vectors 
forms a basis for the points generated by the linear con- 
gruential generator. The objective of the lattice test is 
to start with an initial basis and reduce it through elemen- 
tary row operations to a so-called optimal basis. The 
characterization of the sequence is given by the ratio of 
the longest side to the shortest side of the optimal hyper- 
cube defined by the reduced basis vectors. 

Ideally, all of the points in n-space will be generated 
so that the optimal lattice basis will have all sides of 
length one and a unit cell volume of one. Since a linear 
congruential generator produces only p points in n-space 
where p is the modulus of the generator, the unit cell volume 
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will be p n For a good generator, the ratio of longest 

to shortest side will be close to one, indicating a sem- 
blance of regularity of distribution of the points. An 
ad hoc rule is to summarily reject any generator whose side 
ratio is greater than two . 

A natural starting basis for the case of n = 2 is given 
by the vectors 



(l,a) and (0,p) 

where p is the modulus of the generator. These basis vectors 
are derived as follows: 



x i+i ■ a x i • kp 



Vi " a x j ' lp 



for some integers k and 1, so 



(Xj,Xj +1 ) - ( x j/ x i+1 ) = (L,aL-(l-k) p) 

where L is the integral distance between and . In 
basis vector form (X^,Xj +1 ) - (X^,X^ +1 ) = L(l,a) - (1-k) (0,p) 
hence the vector difference between any two generated points 
is the sum of integral multiples of the two basis vectors 
(1, a) and (0,p) . This form of a starting basis can naturally 
be generalized to an n-dimensional basis. 
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Marsaglia formalizes this basis concept with a theorem 
for the lattice in 3-space. 

Theorem 8 (Marsaglia [Ref. 91) : Let a. , a„, ..., a be 

12 m 

the set of points in 3-space of the form (x. , x. Ll , x. ,_) 

i i+l i+2 

where the X's are reduced residues of some modulus m gen- 
erated by the congruential generator X i+1 5 a x i + c (mod m) . 
There is no restriction on the integers c, a, or m. If the 
point b = (0, c, ac+c) is subtracted from each of these 
points, then the resulting points, a^-b, a 2 ~b, ..., a m ~b 

all lie on a lattice with unit cell volume m , generated by 
. 2 

the 3 points (1, a, a ), (0, m, 0) and (0, 0, m) . 

The generalization to higher dimensions is straight- 

2 

forward; for example, in 4-space, if b = (0, c, ac+c, a C+ac+c) 
is subtracted from each of the points of the form (x^, x^ + ^, 
x i+2' x i+ 3 ^ ' the resulting a 1 ~b, a 2 ~b, ..., a m ~b all 

3 

lie on a lattice with unit cell colume m generated by the 

2 3 

4 points (1, a, a , a ) , (0, m, 0, 0) , (0, 0, m, 0) and 

(0, 0, 0, m) . 

For comparative purposes, this initial basis is of no use. 
The object of the lattice test is to employ elementary row 
operations on the rows of the initial basis , or unimodular 
transformations to the basis matrix, to acquire an optimal, 
or nearly optimal, basis which can be used to compare 
generators . 

To this end, Marsaglia [Ref. 9] offers the following 
algorithm. 
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Algorithm BEST 2 

Given two points a = (a^a^ . .., a n ) and b = (b x , b 2 , ... , 

b n ) in n-space 

Bl. If b is shorter than a, interchange b and a. 

B2. Replace b by b - La, where L is the integer closest 
to 

ab'/aa' = la^/Ia 2 

B3. If the new b is longer than a, stop. Otherwise, go 
to step Bl. 

BEST2 starts with the initial basis vectors defined 
above. Upon termination, the first row contains the shortest 
reduced basis vector and the last row contains the longest 
reduced basis vector. The ratio of these two basis vectors 
characterizes the lattice produced by the generator in 
question. 

To further understand the implementation of the lattice 
test using BEST2, Marsaglia [Ref. 9] gives the following 
example for the case of n = 4 . 

Algorithm N 

Nl. Start with a basis a^ a 2 , a 3 , a 4 ordered in terms of 

increasing (Euclidean) length |a^| <_ |a 2 | £ |a 3 | <_ | a 4 1 . 
For the congruential generator T(X) = ax + b (mod m) , 
a basis is (1, a, a 2 , a 2 ), (0, m, 0, 0), (0, 0, m, 0), 

(0 , 0 , 0 , m) . 
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N2 . Apply the BEST2 algorithm to , a 2 then to a^, a^ 

then to a^, a 4 then to a 2 , a 3 then to a 2 , a 4 then to 
a^, a 4 . If these 6 applications of BEST2 do not change 
any of a^, a 2 , a^, a 4 then stop. Otherwise order the 
new basis and repeat this step until the 6 applications 
of BEST2 produce no change. 

N3. If |b-jJ < | b 2 | £ | b ^ | £ | h> 4 | is the basis on which 6 
applications of BEST2 have no effect, use the ratio 
r = | b 4 | / | b.^ | to characterize the lattice. A ratio r 
near 1 means a good lattice and a large r means a bad 
lattice. Typically, good lattices have r < 2 while 
r > 3 might be defined as a bad lattice. 

The implementation just presented has become popular for 
implementation on digital computers. Later in this section, 
the inherent computational problems of this algorithm will 
be pointed out. 

In Marsaglia's own terminology the repeated application 
of the BEST2 algorithm results in a nearly optimal basis. 

In certain cases it arrives at the optimal, or fully reduced 
basis. In certain cases it arrives at the optimal, or fully 
reduced basis. Preceeding the work of Marsaglia, Beyer, 

Roof, and Williamson [Ref. 4] formulated an algorithm for 
obtaining the optimal, fully reduced basis. They approached 
the problem of obtaining the basis which would yield shortest 
vector lengths as a problem in the solution of positive 
definite quadratic forms . In his work on the theory of the 
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geometry of numbers , the f amour German mathematician 
Minkowski established criteria for fully reduced basis 
vectors. Ahrens and Dieter [Ref. 2] also use the Minkowski 
theory in developing their algorithm which will now be 
presented. 

A great deal of preliminary definitions and lemmas are 
necessary to adequately establish reduced bases in the 
Minkowski sense. They will not be repeated here but the 
algorithms themselves will be presented to indicate the 
computations necessary for achieving reduced bases. 

Ahrens and Dieter propose two algorithms, one for the 
case of two dimensions and one for the cases 2 < n £ 6. 

It should also be noted that full reduction of the basis 
vectors in the Minkowski sense has only been proven for 
n £ 4. For the cases n = 5 and n = 6 the sufficiency of 
Minkowski's criteria has not been proven and for the cases 
n > 6 no criteria exist. 

For the two dimensional case the following algorithm 
produces a Minkowski reduced basis for the lattice of a 
linear congruential sequence. 

Algorithm MB2 . 

Ml. Set i = 1, j = 2, and s = 0. 

M2. Set m = [j + (b^b j ) / (b^b^) ] . 

M3. If m ^ 0, go to M5. 

M4. Set s = s+1. If s = 2, go to M7, otherwise go to M6 . 

M5. Set b. = b. - mb. and s = 0. 

D 3 1 
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M6 . Interchange the values of i and j . Go to M2 . 

M7. If |b 1 | > |h> 2 |, interchange and b 2 * 

M8. Return the Minkowski basis b^, b 2 * 

It is assumed that b^ and b 2 are an initial basis as 
in the Marsaglia algorithm. In step M2, the bracket notation 
represents a truncation to the next smallest integer value. 
Products such as (b^b^) represent ordinary vector dot products. 

Full reduction in the cases n = 3 to n = 6 requires an 
auxiliary table, C n , of 5560 coefficients. These may be 
precomputed and stored or they may be computed within the 
algorithm as they are required. The algorithm for these 
cases is as follows. 

Algorithm BMi . 

BM1. Sort the b^ such that |b^| £ | b 2 1 £ | b ^ | £ ... £ | h> n I * 

BM2 . Search for an expression 

m = [i + (b i b.)/(b.b i ) ] 

which is not zero. If no such n exists, set k = 4 
and go to BM4 . 

BM3. Set b. = b. -mb., restore the order and go to BM2 . 

^ ^ 1 n 

BM4 . Set j = n and v = £ c.b. where the c- constitute the 

i=l 1 1 

k-th set of coefficients in the table C n » 

BM5. If Cj = 0 or if the greatest common divisor d = (Cj, 

..., c ) is not 1 set j = j-1 and restart this step, 
n 

(The test of d is easy since d ^ 1 occurs in only four 
possible cases.) If w < bjb j , go to BM7 . 
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BM6 . 



If k points at the last vector in the table C 

n 

(i.e., if k = 26, 80, 402, 4952 for n = 3, 4, 5, 6) 
go to BM8 . Otherwise set k = k+1 and go to BM4 . 

BM7 . Set bj = v restore the order |b^| £ | b 2 | £ ••• £ |b n | 
and go to BM2 . 

BM8 . Return the Minkowski basis b. , b_, ..., b . 

12 n 

The table, C , is used to insure that there remain no 
n 

unimodular transformations which will further reduce a given 
pair of basis vectors. In many cases the basis vectors can- 
not be further reduced. Here the Marsaglia algorithm and 
the Ahrens and Dieter algorithm would stop at identical bases. 
The Ahrens and Dieter algorithm proceeds, however, when 
further reduction is possible. Ahrens and Dieter offer no 
theory indicating which multipliers require the additional 
Minkowski reduction. 

B. COMPUTATIONAL IMPLICATIONS 

Both algorithms just presented are not so complicated 
that a programmer would have difficulty coding them for exe- 
cution on a digital computer. Ahrens and Dieter provide a 
fairly straight forward method for constructing the table of 
coefficients, C n . Marsaglia has provided examples of his 
own algorithm worked out, ostensibly by hand, for some two 
dimensional cases. 

The real problems arise when it is realized that except 
for very small moduli, existing digital computers cannot 
express the values which arise during the course of the 
algorithm's execution. 
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In detailing the algorithms, no mention was made of 
reducing the integer values by the modulus. In fact, in 
the initial basis the modulus appears as the only non-zero 
element in n-1 of the basis vectors. The first vector simi- 
larly contains successively higher powers of the multiplier 
which can be assumed to eventually exceed the modulus in 
magnitude . 

The modulus for a random number generator is typically 
chosen to be a value very close to the word size of the host 
computer. In attempting to implement the lattice test, a 
programmer is immediately confronted with the problem of 
expressing the values which exceed the modulus and conse- 
quently exceed the computer's word size. Knuth [Ref. 8] 
has advised that at least 90-bit integer arithmetic is re- 
quired to accurately compute the desired results for moduli 
35 

of the order 2 . Needless to say, this requirement exceeds 

the capacity of existing general purpose computers. One 
known method to alleviate the problem is to acquire a soft- 
ware package which handles unlimited precision computations. 

The computational problems are certainly a hinderance, 
but not insurmountable. It will now be shown that all of 
these problems are a result of the inappropriate formulation 
of the lattice test. By avoiding the over-sophistication of 
the previous development, useful results can be obtained 
fairly simply. 



37 



C. THE LATTICE TEST IN A FINITE FIELD 

In Chapter II it was shown that by viewing the generated 
points as elements of a finite field, interesting properties 
of the structure of linear congruential sequences could be 
revealed. The use of the finite field construct was not a 
gimmick to achieve the results of that chapter. It is the 
only rational construct to use in examining these finite, 
periodic sequences. 

All of the computational difficulties with the lattice 
test as described above arise from the decision of the authors 
to treat linear congruential sequences as infinite sequences 
of integers, or worse yet, as the field of real numbers. 

The modulus m (composite) or p (prime) is not an element of 
the field and cannot legitimately be expressed in the field. 
Likewise, the intermediate computations which result in 
values greater than the modulus represent quantities not in 
the field which is to be characterized. 

The same concepts of finite field arithmetic which were 
applied to the individual elements of the field in Chapter II 
are readily extensible to the finite field of vectors in a 
finite n-dimensional space. Arithmetic in a finite vector 
space over a field will be examined and an intuitive develop- 
ment of the lattice test will now be presented. 

Conceptually, an n-dimensional space defined on a finite 
field with p elements contains p n distinct points or p n vec- 
tors emanating from the origin (taken to be the vector additive 
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identity (0,0, . .., 0)). Arithmetic defined on these vectors 

possesses the field properties relating to the operations of 
addition and multiplication, that is the properties of clo- 
sure, distributivity , commutativity, and the additive and 
multiplicative identities. 

A linear congruential random number generator with primi- 
tive root multiplier and prime modulus produces a periodic 
sequence of (p-1) integers. Although the element zero does 
not appear explicitly, it will be agreed to include it to 
satisfy the requirement of an additive identity for field 
arithmetic. Hence the p elements (0,1,2, ..., p-1) are 
produced in some permuted order. 

The multidimensional properties of a linear congruential 
sequence are derived from considering consecutive n-tuples 
of the generated elements. Although the n-dimensional finite 
vector space defined on a field of p elements contains p n 
vectors, the set of n-tuples derived from a linear congruen- 
tial sequence contains only p of these vectors, that is, 
the vectors 



({x 0 ,x 1 . 






{x r x 2 . 



,X n ), 



* X p-l' X 0' • ■ • ' X n-2^ 



These p vectors satisfy the additive field property, that is, 
the sum of any two vectors is also a vector in the field. The 
additive identity is defined since for any vector, there is 
another vector such that their sum is (p,p,...,p) = (0,0,..., 0) 
(mod p) . 
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The lattice test need only be concerned with the p vec - 
tors generated by the linear congruential sequence regard- 
less of the dimension of the space being considered. The 
optimal hypercube in n-dimensional space will be defined by 
n vectors of the p possible vectors. The lengths of sides 
of the hypercube will be the difference of the vectors 
defining adjacent corners of the hypercube, and by the pro- 
perties of finite field arithmetic,- this vector difference 
will also be one of the p possible vectors. 

The previous development of the lattice test required a 

basis of n linearly independent vectors to span the space. 

This is a requirement only when all p n vectors are to be 

spanned. Since only p^ possible vectors are to be spanned 

in this formulation, only one basic vector is necessary. 

n _ ^ 

A natural starting basic vector is then (l,a,...,a ; (mod p) . 

It is clear that this vector is one of the p possible vec- 
tors and all other vectors are integral multiples of this 
basic vector. Reduction modulo p is carried out since the 
field operations are defined to be modulo addition and modulo 
multiplication. Problems of overflow are avoided since the 
vector elements never exceed p in magnitude. 

The question of a new lattice test algorithm remains. 

Which integral multiples of the basic vector constitute the 
sides of the smallest n-dimensional hypercube formed by the 
p possible vectors? 

The two dinensional case will be considered. In the 
finite two dimensional space over the prime field of 
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2 

characteristic p, there are p possible vectors. Given a 
primitive root/prime modulus random number generator 
X i+1 = a (mod p) , the full period will consist of a 
permutation of the (p— 1) field elements. Again, it will 
be agreed to include the additive identity element 0. 

Taking successive, overlapping pairs of the elements of the 
linear congruential sequence results in p possible vectors 
in the finite space of p vectors. Since the integer 1 
appears in the sequence, the basic vector will be chosen 
to be (l,a) . The order in which the vectors are generated 
is immaterial. By systematically choosing the vectors, the 
lattice spanned by the p vectors can be constructed without 
lengthy search. 

Assuming the primitive root multiplier to be less than 
^p, the next vector of interest is the vector (2,2a) = 2(1, a). 
Clearly, this vector is one of the p possible vectors and it 
can also be seen that it is oriented in the same direction as 
(l,a) but with twice the magnitude. Another way to view this 
situation is to treat (l,a) and 2(1, a) as points which lie 
on the same line emanating from the origin (0,0), each point 
maintaining a constant distance from the preceeding point. 

To determine how many points lie on this line it is 
necessary to determine which value causes the quantity 
k^a to exceed the modulus p and consequently be reduced 
(mod p) . With k^ reduced (mod p) , the vector k 1 (l,a) (mod p) 
is no longer oriented in the same direction as the vectors 
(1, a) , 2(1, a), ..., (k 1 ~l) (l,a) . The value k 1 is simply 
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p/a where the brackets indicate the greatest integer 
function . 

Beginning with the point k^(l,a) (mod p) another line is 
formed parallel to the first line consisting of the points 

{k 1 (l,a), (k^l) (l,a) , (k 2 *-l)(l,a)} (mod p) . 

The value k^ is 2p/a . This procedure may be continued for 
successive parallel lines until all p points have been 
covered. 

It is, of course, not necessary to generate all of the 
p points nor is it necessary to determine how many parallel 
lines are contained in the space. To find the shortest 
vector emanating from the origin to one of the parallel lines 
requires at most a computation of the first element of the 
line. This bound is actually much greater than the number 
of computations required in practical situations. Assuming 
that the vector to a point k^(l,a) (mod p) has been computed, 
all starting points k^(l,a) (mod p) which lie on the same 
line need not be considered since they are parallel to the 
vector and obviously are greater in magnitude. 

Some examples will now be shown to clarify the procedure 
described above. 

Example 1 . = a (mod p) ; a = 3, p = 17. 

The basic vector is (1,3). k^ = 17/3 = 6, k 2 = 34/3 = 12 

> 3. Only two vectors need be considered. The Euclidean 
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length of the vector (1,3) is (10) 1//2 . The length of the 

1/2 

vector (6,1) is (37) . To complete the search for the 

parallelogram with shortest sides, the vector (6,1) - (1,3) 

= (5,2) is computed and its length is (29) 1//2 . The side 
ratio of the optimal lattice parallelogram is (29) 1//2 / (10) 1//2 
= 1.70. 

Example 2 . X^ + ^ = a (mod p) ? a = 7, p = 17. 

The basic vector is (1,7) . k ^ = 17/7 = 3, k 2 = 34/7 = 5, 

k^ = 51/7 = 8 > 7. The respective lengths are (50) ^ //2 , 

(25) and (26) l^ 2 . The vector (3,4) with the smallest 

1/2 

length (25) will be one side of the optimum lattice 

parallelogram. The distance from (3,4) to (5,1) , the next 

smallest, will be computed to determine if it will yield a 

1/2 

smaller side. Since this length is (13) , the optimum 

lattice parallelogram has side ratio (25) ^ 2 / (13) = 1.39. 

Example 3 . X^ + ^ = a (mod p) ; a = 7^ = 16,807; p = 2 2 ^ - 1. 

The basic vector is (1,16807)'. k^ = 2 2 ^ - 1/16807 = 127774 

> 16807. The respective lengths are (282475250) ^ 2 and 
(16521383917) 1/2 . The length of the vector (127774,13971) 
to (1,16807) is (16333982425) 1</2 . The side ratio of the two 
shortest sides is then (16333982425) l/ 2 / (282475250 ) ^^ 2 = 4.74. 

Two interesting results should be noted here. First, 
this generator has less than optimal two dimensional lattice 
characteristics and second, the previously reported results 
for this generator used the ratio 
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(16521383917) (282475250) = 7.60 from Marsaglia's 
algorithm. This is evidence of the Marsaglia algorithm 
BEST2 not achieving an optimal basis. 

This chapter has unfortunately stopped tantalizingly 
short of producing an algorithm for use in the many cases 
of interest. The last chapter will outline the course of 
future work in this direction. 
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V. THE SPECTRAL TEST 



The spectral test for linear congruential random number 
generators is directly related to the lattice test. The 
development of the spectral test comes from an entirely 
different theoretical basis, however. 

The application of the spectral test to linear congruen- 
tial random number generators began with the work of Coveyou 
and Macpherson. Knuth [Ref. 8] published a computational 
algorithm which, as stated in the previous chapter, requires 
90-bit integer arithmetic to achieve accurate results. A 
brief historical outline of the basis and development of the 
spectral test will be given. The test will then be compared 
and contrasted with the lattice test. 

A. HISTORICAL DEVELOPMENT 

The theoretical basis of the spectral test for linear 
congruential sequences is a theorem by H. Weyl in 1916. The 
theorem was concerned with the distribution of sequences of 
deterministic numbers reduced modulo 1. Although the sequences 
to be examined were deterministic in nature, they were also 
assumed to be infinite in length and of infinite precision, 
in representation. 

Assuming an infinite sequence of real numbers, say Xq,X^, 
..., X , ..., on the half-open interval [0,1), select a number 
Y also on the interval [0,1). Let N y be the cardinal number 
of the first N values of {x} which are contained in the 
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subinterval [0,Y). The sequence {X} is then said to be 

equidistributed modulo 1 if lim N /N = Y for all values of 

N-*-eo y 

Y. Weyl's theorem extends this definition to an arbitrary 
number of dimensions . 

A 1963 paper by Franklin [Ref. 6] applied Weyl's theorem 
to sequences of pseudo-random numbers. Franklin was con- 
cerned with the equidistribution of the deterministic se- 
quence X^ = 0 1 (mod 1) for i = 1,2,.... He proved that this 
sequence is completely equidistributed for almost all 0 > 1 
and further proved that 0 must be a transcendental number. 
Franklin's work was also developed on the basis of infinite 
precision of the values X^ (mod 1) . 

Although they do not cite Franklin's work as a precedent, 
Ahrens and Dieter [Ref. 3] have proposed a FORTRAN implemen- 
tation of a random number generator using real arithmetic 
and the golden section number as a multiplier. This number 
is, of course, irrational but not transcendental. They 
claim the equidistribution properties of this generator are 
derived from the irrationality of the multiplier. 

The Coveyou and Macpherson development was the first 
time that Weyl's theorem was applied to finite, periodic 
sequences of linear congruential generators. Knuth amplified 
the development of the Coveyou and Macpherson test and 
offered a complete algorithm for its implementation. 
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B. THE ALGORITHM 



Weyl's theorem provides necessary and sufficient condi- 
tions for equidistribution of infinite sequences of deter- 
ministic numbers reduced modulo 1. To introduce Weyl's 
theorem, let X^, i = 1,2,... be a sequence of n-tuples, that 
is, X. = (X . , X.. ., ..., X . .) where the X. . all lie 

X U / X X f X Xl*" X t X J / X 

in the interval [0,1) . The X^ are n-dimensional vectors in 
the n-dimensional unit hypercube in Euclidean space. 

Weyl's Theorem .^ The sequence of points X^,^/ ••• is 
equidistributed modulo 1 if and only if 



lim 4 

p-*c 



p ^ ex P C"2Tri , j 



Vl,j + •••- = 



for all vectors (q^ ,q^, . . . ,q n _^) with integer elements not 
all zero. 

The proof of this theorem is rather lengthy and will not 
be repeated here. Coveyou and Macpherson [Ref. 5] applied 
this theorem to the case of linear congruential random number 
generators and provided computational details for an algorithm. 
Knuth's development [Ref. 8] of the algorithm will be 
outlined. 

To begin, the spectral test is only applicable to full 
period generators, that is, properly formulated mixed 



"^Jansson, Birger, Random Number Generators , p. 156, 
Amlquist and Wiksell, Stockho lm, 1966 . 
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congruential generators with composite modulus or primitive 
root/prime modulus generators. Although Knuth treated the 
mixed congruential case, the primitive root/prime modulus 
case will be substituted in this brief outline. 

Let the random number generator be X^ + ^ 5 a (mod p) . 
The finite Fourier transform of this sequence is then 



P-1 



-2iTi 



f (s, ,s 0 , . . . ,s ) = - E exp ( 

12 n P k=0 P 



(s l X k +s 2 X k+l + * • ,+s n X k+n-l ) ) 



p k=0 



exp ( 



-2iTi 

P 



(s (a) X k ) ) 



where 



/x , n-1 

s(a) = s^ + s 2 a + " ' ' s n a 

Since all values of k appear in the sequence and their order 
is immaterial, the following may be substituted: 

1 P”1 -2ui 

f (s. ,s_ , . . . ,s ) = — I exp ( — - — s (a) k) . 

1 z n P k =o p 

Noting that this represents the sum of a geometric series , 
the basis formula is 

f (s^, S 2 / • • . / s n ) = 6(s(a)/p) (1) 

where 6 (X) = 1 if X is integer and 0 otherwise. 



48 



Application of discrete probability theory yields the 
result that the joint probability of any n-tuple (x k , X k+1 , 
^k+n-l^ hence the value f (s^, s^ , . . . , s n ) /p 

may be interpreted physically as the amplitude of the 
n-dimensional complex plane wave 



«(ti, t 2 , . . . ,t ft ) = exp ( 



27ri 



( S . t . t . . • t s t ) ) • 
11 n n 



A wave number may assigned to this plane wave corresponding 
to the frequency of the wave. The wave number is 



v 




.+s 



n 



2 } l/2 



for s k p/2. 



In a truly random sequence no waves should be present 

except a constant wave of frequence zero. Any generator 

which produces a swquence whose transform yields nonzero 

frequencies can be considered to be nonrandom. To summarize 

2 

the spectral test Knuth states: 

"If v is the smallest nonzero value of the 
n 

wave number ... for which f (s^ , s 2 , . . . , s n ) / 0 
in a linear congruential sequence with maximum 
period, then the sequence X^/m, x^/m, X 2 / 111 , . . . 
represents a sequence of random numbers uniformly 
distributed between 0 and 1, having 'accuracy' 



9 

Knuth, D.E., The Art of Computer Programming, Volume 2: 
Seminumerical Algorithms , p. 85, Addison-Wesley , N. Y . , 1969 . 
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or 'truncation error' 1/v with respect to 
the independence of n consecutive values of 
the sequence averaged over the entire period." 

Equation (1) represents the spectrum of the linear con- 
gruential sequence. It can be seen from (1) that 
f (s^, d 2 , . . . , s^) = 0 except when 

s (a) = s 1 + s 2 a + ...+ s n a = 0 (mod p) 

in this case f (s^ s 2 , . . . ,s ) = 1. Again directly quoting 
Knuth 2 : 

"... for linear congruential sequences of 
maximum period, the smallest nonzero wave 
number in the spectrum is given by 

v = min (s 2 + s 2 + ...+ s 2 ) 1/2 
n 1 z n 

where the minimum is taken over all n-tuples of 
integers (s^, s 2 , • • • , s n ) satisfying (2)." 

Knuth then proceeds to develop an elaborate computational 
scheme to implement the spectral test. As with the lattice 
test algorithm of Ahrens and Dieter, the test algorithm is 
formulated in terms of the solution to a positive definite 
quadratic form. To find the minimum value of the wave 



2 Ibid. , p . 85. 
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number ( s ^ + s 2 + ...+ such that 

n “ 1 

s(a) = + S 2 & + ... + s n a = 0 (mod p) the problem may 

be reformulated' as : 

find the minimum value of the quantity 



(a ll S l +a 12 S 2 + * * * +a ln S n ) + 



* + (a nl s l +a n2 s 2 + * * * +a nn S n ) 



This results in the positive definite quadratic form 



( S ^ , S 2 t • • • / S^) ' ^ ( S l /S 2^ * * • /S n ) 

where (s^, S 2 / • • . / s^) is a column vector, not all zero, and 
A is any nonsingular matrix of coefficients. 

C. THE SPECTRAL TEST IN A FINITE FIELD 

The computational algorithm for the spectral test is 
somewhat more involved than the algorithm for the lattice 
test. Knuth's algorithm, however, is very carefully written 
and can very easily be implemented assuming the availability 
of multiple-precision arithmetic subroutines. 

In view of the finite field approach of this thesis, the 
basic assumptions of the spectral test appear entirely in- 
appropriate. To begin, consider the discrete probabilistic 
basis of the spectral test. In one dimension, the probability 
measure assigned to each element is correctly 1/p. However, 
in higher dimensions, the probability measure for each n-tuple 
in the space is taken to be l/p n . Since only p distinct 
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* 



n- tuples are generated by a linear congruential generator, 
the correct measure should be still 1/p. This point is more 
clearly seen when it is recalled that all n-tuples generated 
by a primitive root/prime modulus generator in finite dimen- 
sional space are integral multiples of the basic vector 
(l,a,a 2 , . . . ,a n ^) . That is, the set of all possible n-tuples 
is 



{k(l,a,a 2 , . . . ,a n 1 ) (mod p)}, k=l,2,...,p-l 

Similarly, the theoretical basis of the spectral test is 
the finite Fourier transform defined only in the complex 
field. When viewing the sequence generated by the primitive 
root/prime modulus random number generator as a finite field, 
a more appropriate transform is available, specifically the 
number theoretic transforms as developed by Rader [Ref. 12], 
Agarwal and Burris [Ref. 1] , and Pollard [Ref. 11] . The 
finite Fourier transform has been used since in the complex 
field a primitive root of unity exists, namely exp(2"rrik). 

In finite fields, any primitive root of the prime modulus is 
a primitive root of unity and consequently a transform (gen- 
erating function) is defined. For the case of primitive 
root/prime modulus generators, the spectral test may be more 
appropriately defined in terms of these transforms. 

As in the lattice test, only p possible n-tuples can 

be generated, hence the search for the minimum value of the 

22 2 1/2 

wave number (s^ + s 2 + • • ♦ + s n ) such that 



52 



s (a) - s ]_ + s 2 a + ••• + s^a 11 ^ = 0 (mod p) becomes impossible. 
The following theorem demonstrates this fact. 

Theorem 9 : For the primitive root/prime modulus random 

number generator X i+1 E a X i (mod p) no solution, to the 
basic congruence 

n— 1 

s(a) = s^ + S 2 a + ... + s n a = 0 (mod p) 

exists in the finite field of characteristic p. 

Proof : Referring to the development of finite fields in 

Chapter II and Chapter IV, Section C, the only possible 

wave in the spectrum of the primitive root/prime modulus 

random generator must be of the form k(l,a,...,a n (mod p) 

for some integer k such that 0 < k <_ p-1. Since p is prime, 

n _ 

it has no integral divisors, hence (l+a+...+a x ) must be 
congruent to 0 modulo p. This cannot happen for any n such 
that 0 < n £ p-1. Let n = 1, then 

(1+a) = 0 (mod p) 

since this implies that a = p-1 or a = -1 in primitive mark 
notation and ±1 are not, by definition, primitive roots of 
any prime p > 2. Let n = p-1, then 
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( 1 + a + . . . + ^ ) 



P-1 . 

- Z a 1 = 1 - a^/l-a 
i=0 

= 1 - a/l-a 



=1^0 (mod p) . 



Let 1 < n < p-1, then 

n-1 

(1 + a + . . .+a n ) = = a 1 = l-a n /l-a % 0 (mod p) 

i=0 

since for this to be congruent to 0 would imply that a 11 = 1 
which contradicts the primitivity of a for which p-1 is the 
smallest positive integer such that a^ 3 ^ = 1 (mod p) . Q.E.D. 

The spectral test, when considered as a special case of 
the lattice test, is still useful when characterizing linear 
congruential sequences. Unfortunately, the formulation of 
the test is, at best, shakey. In simpler terms, the spec- 
tral test is designed to determine the minimum Euclidean 
distance between consecutive hyperplanes containing points 
of the generated sequence. In contrast to the lattice test, 
only the shortest side of the hypercube is of interest, rather 
than the ratio of longest to shortest side. The spectral 
test may be most simply restated for primitive root/prime 
modulus random number generatros as : 
find the minimum value of 
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• • • 



) (mod p) 



v 2 = k ( l 2 + a 2 + 
n 



+ a 



2 (n-1) 



over all k such that 0 < k <_ p-1. 

As in Chapter IV on the lattice test, this chapter will 
stop painfully short of producing a general purpose algorithm 
for a new spectral test. As the conclusions will reiterate, 
future work will most certainly produce such an algorithm. 
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VI. CONCLUSIONS 



The development of primitive root/prime modulus sequences 
as finite fields has opened an entirely new perspective on 
the characterization of congruential random number generators. 
The results established in the proceeding chapters have been 
confined to the class of primitive root/prime modulus genera- 
tors but it is felt that the results are readily extendible 
to mixed generators with composite moduli. 

The lattice test and spectral test have been shown to be 
equivalent in their power to characterize sequences, the 
choice of one test over the other is merely a matter of 
taste. While the theoretical, geometric, and intuitive 
appeal of both tests are appealing, their implementation has 
suffered from over-sophistication. The sparse results pre- 
sented here were worked out by hand calculator. Several 
other generators have also been examined by hand and the 
results agree with published lattice test results. Some 
cases of published results did not work out as simply since 
the multipliers were greater than the square root of the 
modulus and the multiplicative inverses, which would have 
produced the correct results, were too difficult to find. 

Further development of the results presented here will 
lead to criteria for selecting multipliers which will produce 
nearly optimal properties of the sequence as characterized 
by the lattice and spectral tests. 
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The lattice and spectral test algorithms can be rede- 
veloped completely based on finite field properties. The 
implementation of such algorithms will be quite natural for 
digital computers since all arithmetic on digital computers 
is necessarily confined to finite field arithmetic due to 
the finite capacity of the registers. ("Real arithmetic", 
or floating point, is confined to a finite set of the 
rationals . ) 

The implications of the results presented here are far- 
reaching. On-going research should bring them to fruition. 
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